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ABSTRACT 

We present direct numerical simulations of the equations of compressible magnetohydrodynamics 
in a wedge-shaped spherical shell, without shear, but with random helical forcing which has negative 
(positive) helicity in the northern (southern) hemisphere. We find a large-scale magnetic field that 
is nearly uniform in the azimuthal direction and approximately antisymmetric about the equator. 
Furthermore, the large-scale field in each hemisphere oscillates on nearly dynamical time scales with 
reversals of polarity and equatorward migration. Corresponding mean-field models also show similar 
migratory oscillations with a frequency that is nearly independent of the magnetic Reynolds number. 
This mechanism may be relevant for understanding equatorward migration seen in the solar dynamo. 
Subject headings: MHD - Turbulence 
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1. INTRODUCTION 

Large-scale magnetic fields with fascinating quasi- 
regular spatio-temporal behavior are ubiquitous in solar 
and stellar settings. Understanding the mechanisms for 
the generation of such fields and their spatio-temporal 
variations is still a major challenge for dynamo theory. 
The solar magnetic field has three particularly important 
features: quasi-regular oscillations, reversal of polarity 
and equatorward migration. Direct numerical simula- 
tions of solar-like convective dynamo s have been able to 
generate large-sca le magnetic fields ([Brown et al.l 120071 : 
iBrown et al.lboiOl) which in some cases show oscillatory 
behavior, but the fields exhibit eith er rather weak equa- 
torward migration at high latitudes dGhizaru et alJlToiOf) 
or anti-solar (i.e. poleward) migration ( Gilmanl 119831 : 
IKapvla et a!]l20Toft . 

A useful tool for studying these dynamical phe- 
nomena is mean-field (MF) electrodynamics (e.g. 
Krause fe Radlerl 119801 : IBrandenburg fc Subramanianl 
200.^ . where the effects of turbulence are character- 



ized by turbulent magnetic diffusion and an a effect. 
According to MF theory, equatorward migration is ex- 
pected if there is negative radial shear accompanied by 
a positive (nega tive) a effect in the n orthern (south- 
ern) hemisphere (|Krause & Radlcr 1980). Direct numer- 
ical simulations (DNS) of helical turbulence with shear 
have confirmed the presence of migratory dynamo waves 
(|Brandenburg et alJ200ll IKapvla fc Brandenburjl2009h . 
It is, however, unclear whether this is really what is 
going on in the Sun, since there the layer with nega- 
tive radial shear is rather thin and on l y con centrated 
near the surface (see, e.g., IBrandenburg! (|2005j ) and ref- 
erences therein). The other alternative is that merid- 
ional circulation might change the direction of migration 
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( Choudhu ri et al.lll995l) , but evidence for this has not yet 
been seen in DNS. 

In this Letter we present a completely different mech- 
anism for polarity reversal and equatorward migration 
of dynamo activity. In the context of MF models this 
mechanism is connected with th e antisymmetry of the 
profiles of a across the equat or (|Rudiger fc Hollerbachl 
12004 IBrandenburg et aLll2009D . We demonstrate the op- 
eration of this mechanism in DNS of the equations of 
magnetohydrodynamics (MHD). Our model consists of 
a spherical wedge-shaped shell in which the turbulence 
in the fluid is maintained by a random helical forcing. 
Motivated by the Sun, we choose our forcing to have op- 
posite signs of helicity in the two hemispheres (negative 
in the north and positive in the south). We emphasize 
that, even though our model does not explicitly include 
convection, stratification and rotation, the helical forcing 
used here does partially model these features implicitly. 

Our model shows large-scale magnetic fields in excess 
of the equipartition value. More importantly, we find 
oscillations of the magnetic field which show opposite 
signs in different hemispheres with periodic reversals of 
polarity. Furthermore, the magnetic field develops at 
higher latitudes and migrates equatorward where the two 
different polarities of magnetic field annihilate and the 
cycle repeats itself as shown in the top panel of Fig. [TJ 
To our knowledge, such dynamical features of the large- 
scale magnetic field have not been observed earlier in 
DNS of MHD turbulence. Below we introduce our model, 
discuss its oscillatory solutions, and briefly compare our 
DNS results with those obtained from corresponding MF 
models. 

2. THE MODEL 

In our simulations, we solve the equations for com- 
pressible MHD in terms of the velocity U, the logarith- 
mic density lnp, and the magnetic vector potential A, 



D t U = — Cg V lnp + - J x B 

P 

D t lnp=-V-Z7, 

d t A = U x B + t]V 2 A, 



(1) 

(2) 
(3) 



B t /B 




l/T 



Fig. 1. — Space- time diagrams of the azimuthal component of the 
large-scale magnetic field for DNS over both the hemispheres (top 
panel), DNS over only the northern hemisphere with an antisym- 
metric condition at the equator (middle panel), and MF simulation 
over only the northern hemisphere with an antisymmetry condition 
at the equator (bottom panel). 
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Fig. 2. — Orthographic projection of the toroidal magnetic field 
£?0 at r = 0.85 for Run S5. The projection is tilted by 15 degrees 
towards the viewer. 

where F visc = (p,/ p)(V 2 U + ^ VV • U) is the viscous 
force, ji is the dynamic viscosity, B = V x A is the 
magnetic field, J — V x B/(M) is the current density, 
po is the vacuum permeability, c s is the (constant) speed 
of sound in the medium, r\ is the magnetic diffusivity, 
and D t = d t + U ■ V is the advective derivative. Our 
computational domain is a spherical wedge with radius 
r G [7*1,1*2] symmetric about the equator with colatitude 
9 £ [O, 7r — O] and azimuth <f> G [0, $]. The radial, merid- 
ional and azimuthal extents of our domain are respec- 
tively, L T = r 2 — tt, Lq = r 2 (7r — 20) and = r 2 <I>. All 
lengths are measured in the units of r 2 . Our main results 
do not depend on our choices of and $. 
In Equation [T] f(x,t) is an external white- 



in-time random helical forcing constructed 
using the Chandra se khar- Kendall functions 
(|Chandrasekhar fc Kendal] I1957D as described be- 
low. In spherical coordinates a helical vector function 
can be expressed in terms of a scalar potential function: 

1>(P(t),l(t)Mt)) =Rcz e (Pr)Y e m (eA)eM^rn(t)l (4) 

with zt(/3r) = agji{j3r) + bgni{j3r). Here jg and ng are 
spherical Bessel functions of the first and second kind 
respectively, ag and bg are constants determined by the 
boundary conditions and £ m is a random angle uniformly 
distributed between and 2ir. The helical forcing / is 
then given by the equation Vx/ = (3f, where / = T+S, 
T = V x (etp), and S = fi^V x T, where e(t) is a unit 
vector chosen randomly on the unit sphere. As to the 
choice of boundary conditions, wc demand that / is zero 
at the two radial boundaries r — r% and r — r 2 which 
yields the following transcendental equation relating ag, 
bg and j3: 

adtifin) + bingifirx) = atjg{(3r 2 ) + b e n e (/3r 2 ) = 0. (5) 

We first construct a table of values of m, £ and (3 in 
the following way. As we use periodic boundary condi- 
tions along the azimuthal direction, m m i n = 27r/$. We 
choose m = pm m ; n , and for a fixed m we choose i to 
be odd, £ = 2(m + q) + 1, because we want the forcing 
to go to zero at the equator. Here p and q are integers 
which range between 3 to 5. For a fixed £ and m we 
solve Eq. ([5]) by Newton-Raphson method and list the 
solutions which have 3 to 5 zeros in the domain. To ran- 
domize the resulting forcing we randomly choose a triplet 
of m, £ and j3 from the table. We also randomize the unit 
vector e on the unit sphere. Two different signs of he- 
licity are imposed by choosing negative (positive) j3 in 
the northern (southern) hemisphere. The choice of pa- 
rameters implies that we have a scale separation between 
3 to 5 in our simulations. Our results are fairly robust 
under the change of different parameters of forcing. We 
need scal e separation of 3 or more to excite a lar ge-scale 
dynamo (jHaugen. B randenburg, "fc Doblerl 120041 ) . which 
invariably shows oscillations and equatorward migration. 
Ou r simulations are pe rformed using the Pencil Code0; 
see iMitra et"aL1 (|2009fl for details regarding the imple- 
mentation of spherical polar coordinates. 

We use periodic boundary conditions along the az- 
imuthal direction and set the normal component of the 
magnetic field to zero on all other boundaries (perfect- 
conductor boundary condition). This is implemented by 
setting the two tangential components of A to zero. As 
an estimate of the characteristic Fourier mode of the forc- 
ing we define kf = w rms /u rms , (column 8 of Table [1} 
where w rms and w r ms are the rms values of small-scale 
vorticity and velocity, respectively. We introduce the 
fluid and magnetic Reynolds numbers as Re = u TTas /vki 
and Rcm = Ui ms /vkf> respectively. Here, v = fi/po is 
the mean kinematic viscosity, where po is the initial and 
the mean density in the volume. A representative list of 
parameters is given in Table [T] 

3. RESULTS 

We start our simulations with a random seed magnetic 
field of no particular parity about the equator. After a 

5 http : //pencil-code . googlecode . com 
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Run 


Grid 


Le 


L$ 


-Srms/-Scq 


Re 


Re M 


fc f /fci 


x 10 5 


tj x 10 5 


?7t x 10 5 


lj c x 10 3 


T x 10" 3 


^max 


SI 


32 x 64 x 32 


7r/5 


tt/10 


0.88 


5 


12 


3 


5 


2 


5.3 


2.5 


0.18 


~ 10T 


S2 


64 x 128 x 64 


tt/5 


tt/10 


0.79 


8 


21 


4 


3 


1.2 


6.2 


2.5 


0.16 


~ 5T 


S3 


32 x 64 x 64 


tt/5 


tt/5 


1.16 


2 


4 


7 


5 


2 


3.6 


2 


0.27 


~ 10T 


S4 


64 x 128 x 128 


tt/5 


vr/5 


1.04 


4 


10 


7 


2 


1 


4.9 


2.6 


0.2 


~ 5T 


S5 


64 x 128 x 128 


97T/10 


tt/2 


2.03 


7 


13 


7 


2 


1 


4.4 






~ T. 



TABLE 1 Summary of our parameters including grid size, the meridional and azimuthal extents of our domain, rms value of the 
azimuthally averaged field, Reynolds number Re and magnetic Reynolds number ReM- The forcing amplitude, /amp = 0.2, is chosen 
such that the Mach number is of the order of 0.1, making the flow essentially incompressible. t max the duration of each run. The run S5 

has not been ran long enough to accurately measure ui c . 




Fig. 3. — The frequency of the oscillations multiplied by the 
turbulent diffusion time T = (r/tfcg) -1 as a function of ReM for 
Lg = tt/5. Data from DNS (closed circles) and MF models (open 
circles) are shown. For the MF runs r\t = 1, for the DNS runs rjt 
is given in Table Q] The inset shows uj c X 10 3 (column 12 Table [TJ 
as a function of ReM for the DNS data. 

transient time, of about one turbulent diffusion time, we 
find that a large-scale magnetic field is generated with 
energy of the order of or exceeding the equipartition 
strength in all runs. The magnetic field encompasses the 
whole azimuthal extent of the domain. A contour plot of 
the toroidal component of the magnetic field on a surface 
with constant radius from Run S5 is shown in Fig. [5] We 
define the large-scale magnetic field via averages over the 
azimuthal and radial directions, i.e., B — {B) ri j,, such 
that the resultant magnetic field is solely a function of 
latitude and time. We normalize the magnetic field with 
the equipartition field strength, B cq = (fiopu 2 ) 1 / 2 , where 

u = U — U is the small-scale velocity. The field first de- 
velops at higher latitudes and then with time migrates 
equatorward. In each hemisphere the field shows oscil- 
lations and reversals of polarity. These features can be 
seen in the space-time diagram shown in the top panel 
of Fig. [T] where we plot B^ as a function of latitude and 
time. The principal frequency of oscillations, uj c (col- 
umn 12 of Table [Hand the inset of Fig. [3]) is obtained by 
Fourier transforming the time series of 8^(8, t) in time at 
a given 8 (8 = 7r/20 say) and determining the frequency 
corresponding to the dominant mode. Normalized en- 
ergy in the large-scale magnetic field also shows oscil- 
lations as a function of time, but with frequency 2uj c . 
A characteristic dynamical time scale is the turbulent 
diffusion time corresponding to the length scale Lg de- 
fined by T = (fcg?7t) 1 where kg — 2tt /Lg and for rj% we 
take the expression from the first-order-smoothing ap- 
proximation, rjt = M lms /3fcf. In all our runs we find the 
product uj c T to be of order unity (see the bottom panel 
of Fig.©. 

We note that the radial and azimuthal components of 



the large-scale magnetic field are almost antisymmetric 
about the equator. This allows a further simplification 
of our model by simulating only one half of the domain 
(e.g., the northern hemisphere), while keeping exactly 
the same forcing function (e.g., a forcing that is ran- 
dom, and negatively helical in the northern hemisphere 
going smoothly to zero at the equator), but choosing the 
boundary condition B r = B^ = at the equator. Such 
simulations produce exactly the same oscillations (as can 
be seen by comparing the top and the middle panels of 
Fig. [I) as those obtained in the DNS with both hemi- 
spheres. This implies that these oscillations can be stud- 
ied with half the number of grid points and appropriately 
chosen boundary conditions at the equator. 

Given the large values of the magnetic Reynolds num- 
ber in solar/stellar settings, an important question is how 
the frequency uj c scales with ReM? This question cannot 
be answered from DNS because the magnetic Reynolds 
numbers reached are far from the asymptotic limit of 
large ReM (see Column 7 of Table [lj A way forward 
is to use analogous MF models. The appropriate set- 
ting would be that of an a 2 dynamo with d ynamical a 
quenching (jBlackman fc Brandenburgl I2002T ). which in- 
corporates conservation of magnetic helicity, given by the 
equations 



d t B = V x ( a B) + {i 1 + i lt )V 2 B, 



d t a — —2r/k 



2 ( ct(B ) — rjt(B ■ J) a - a K \ 
f ^ B 2 q Rc M /3 ) 



(6) 
(7) 



where Rcm/3 = and B oq is the equipartition field 
strength. We use fcf/fci = 6 in our MF simulations. In 
view of the discussion above, we solve the MF equations 
in only the northern hemisphere with appropriate bound- 
ary condition at the equator. In the MF approach the 
helical nature of turbulence is modelled by the a coeffi- 
cient («k)- We choose q;k = g{S)ao and r] t = 1. The 
profile function g takes positive (negative) values in the 
northern (southern) hemisphere, going smoothly to zero 
at the equator. This reflects the fact that, according to 
MF theory, the kinetic a effect usually has the opposite 
sign to the mean kinetic helicity. We have used three 
different functional forms for g, namely g = 8 — tt/2, 
g = sin((? — tt/2) and g — tanh(# — tt/2), without any sig- 
nificant change in our results. We need ao > 16 to excite 
a dynamo but once excited the oscillatory and migratory 
properties of the dynamo do not depend on ao. We use 
perfect-conductor boundary conditions along the radial 
direction and our magnetic Reynolds number (changed 
by varying 77) ranges between 10 < ReM < 10 s . We have 
also used domains with larger latitudinal extents than 
those used in our DNS. 
Here we briefly mention a few important outcomes of 
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our MF results relevant to our discussion above: (a) Our 
DNS results - namely oscillatory behavior as well as mi- 
gration towards the equator - are qualitatively repro- 
duced by the MF solutions in the range of parameters 
reported here. An example of the space-time diagram 
produced by our MF runs is shown in the bottom panel of 
Fig- [2 (b) The frequency of oscillations remains almost 
constant as a function of ReM, see Fig. [3l The depen- 
dence of the oscillation period on ReM seen in DNS may 
be related to ReM-dependence of the turbulent magnetic 
diffusivity. A simi lar behavior has been observe d earlier 
in Cartesian DNS (jKapvla fc Brandenbur dl2009T >. (c) To 
show the robustness of our results with respect to the size 
of the domain, we also studied domain sizes extended in 
the meridional direction from Lq = tt/5 to (178/180)tt 
which corresponds respectively to = 72 degrees and 1 
degree. We find that the oscillations and the migratory 
behavior do not change, (d) For the MF model con- 
sidered here the mean value of the large-scale magnetic 
field decreases as Re^ 1 , i.e., the field is catastrophically 
quenched for large values of ReM- In the DNS, however, 
such quenching could be alle viated by magnetic helic- 
ity fluxes across the equator ([Brandenburg et all 120091 : 
iMitra et alJl20ToT ). 

To test the robustness of our simulations with respect 
to the choice of boundary condition in the radial direction 
we repeated our simulations with the vertical field bound- 
ary condition, which makes the two tangential compo- 
nents of the magnetic field vanish at the radially outward 
boundary. These simulations also show oscillations and 
equatorward migration of magnetic activity, but in this 
case the oscillations are less regular and the frequency is 
marginally higher. 

4. CONCLUSIONS 

We have found large-scale fields, oscillations on dynam- 
ical time scales and polarity reversals with equatorward 
migration of magnetic activity in direct numerical sim- 
ulations of helically forced MHD equations in spherical 
wedge domains. Despite its simplicity, it is quite striking 
how our model can reproduce these important features 
of the solar dynamo. As far as we are aware, these fea- 
tures have not been observed earlier in DNS. We have 
elucidated our DNS results by considering analogous a 2 
MF models which support our conclusions. We have fur- 
ther used these MF models to explore magnetic Reynolds 
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numbers that are at present inaccessible to DNS. This 
has enabled us to show that the frequency of the oscilla- 
tions is almost independent of ReM for large ReM- Such 
MF models have been known to have oscillatory solu- 
tion s if a changes sign in the comput a tional domain (see, 
e.g. IBarvshnikova fe Shukurovlll98'7t IStefani fe Gerbeth l 
1 20031; iRudiger fe Hollerbachl 120041: iBrandenburg et al.1 
2009), but their migratory property had not been stud- 
ied before. Antisymmetry of a with depth also pro- 
duces oscillatory solutions rtBaryshnikova fc Shukurovl 
ll987tlStefani fc Gerbethll2003D . but not the equatorward 
migration. 

The helical forcing used in our DNS, with its differ- 
ent signs of helicity in different hemispheres, implicitly 
models only the helical aspect of the effects of rota- 
tion and stratification present in the Sun. Physically, 
a more complete picture should emerge from DNS of 
convect i ye turbulent dynamo as done, for example, b y 
Gilmanl (fl983l) ; iBrown et al.l (120071) :lBrown et all (12010ft: 
Kapyla, et al.l (|2010D ; iGhizaru et al.l (|2010D . Such simu- 
lations also generate differential rotation and lead to an- 
other dynamo mode of operation - the afi dynamo - 
which also produces oscillatory behavior. However, in 
order to get equatorward migration radial shear must be 
negative. Helioseismology has shown that negative radial 
shear exists only near the surface of the convection zone. 
This feature has so far not been reproduced by global 
DNS. Whether or not an a£L dynamo is the dominant 
mechanism operating in the Sun remains unclear. It is 
therefore important to keep in mind that there exists al- 
ternative mechanisms for producing oscillatory behavior 
with equatorward migration, such as the one discussed 
here. 
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